import numpy as np
import pandas as pd
import matplotlib as mpl
import matplotlib.pyplot as plt

from sklearn.tree import DecisionTreeRegressor
from sklearn.ensemble import BaggingRegressor

mpl.rcParams['font.sans-serif']=[u'simHei']
mpl.rcParams['axes.unicode_minus']=False

## 文件存储路径
path1='../data/C0904.csv'
path2='../data/C0911.csv'
filename='X12CO2'

### 原始数据读取
data=pd.read_csv(path1,header=0)
X=data[filename].values

plt.figure(figsize=(10,6), facecolor='w')
plt.subplot(121)
plt.plot(X,'r-',lw=1,label=u'C0904')
plt.title(u'实际数据0904',fontsize=18)
plt.legend(loc='upper right')
plt.xlim(0,80000)
plt.grid(b=True)

data=pd.read_csv(path2,header=0)
X=data[filename].values

plt.subplot(122)
plt.plot(X,'r-',lw=1,label=u'C0911')
plt.title(u'实际数据0911',fontsize=18)
plt.legend(loc='upper right')
plt.xlim(0,80000)
plt.grid(b=True)

plt.tight_layout(2,rect=(0,0,1,0.95))
plt.suptitle(u'如何找到下图中的异常值',fontsize=20)
plt.show()

### 异常数据处理
data=pd.read_csv(path2,header=0)
X=data[filename].values

width=300
delta=5
eps=0.02
N=len(X)
p=[]
# 异常值存储
abnormal=[]
for i in np.arange(0,N-width,delta):
    s=X[i:i+width]
    ## 获取max-min的差值
    min_s=np.min(s)
    ptp=np.ptp(s)
    ptp_min=ptp/min_s
    p.append(ptp_min)
    ## 如果差值大于给定的阈值认为是
    if ptp_min>eps:
        abnormal.append(range(i,i+width))
##获得异常的数据x值    
abnormal=np.array(abnormal).flatten()
abnormal=np.unique(abnormal)
#plt.plot(p, lw=1)
#plt.grid(b=True)
#plt.show()

plt.figure(figsize=(18,7),facecolor='w')
plt.subplot(131)
plt.plot(X,'r-',lw=1,label=u'原始数据')
plt.title(u'实际数据',fontsize=18)
plt.legend(loc='upper right')
plt.xlim(0,80000)
plt.grid(b=True)
 
plt.subplot(132)
t=np.arange(N)
plt.plot(t,X,'r-',lw=1,label=u'原始数据')
plt.plot(abnormal,X[abnormal],'go',markeredgecolor='g',ms=3,label=u'异常值')
plt.legend(loc='upper right')
plt.title(u'异常数据检测',fontsize=18)
plt.xlim(0,80000)
plt.grid(b=True) 

# 预测
plt.subplot(133)
select=np.ones(N, dtype=np.bool)
select[abnormal]=False
t=np.arange(N)
## 决策树
dtr=DecisionTreeRegressor(criterion='mse',max_depth=10)
br=BaggingRegressor(dtr,n_estimators=10,max_samples=0.3)
## 模型训练
br.fit(t[select].reshape(-1,1),X[select])
## 模型预测得出结果
Y=br.predict(np.arange(N).reshape(-1,1))
Y[select]=X[select]
plt.plot(X,'g--',lw=1,label=u'原始值')   # 原始值
plt.plot(Y,'r-',lw=1,label=u'校正值')   # 校正值
plt.legend(loc='upper right')
plt.title(u'异常值校正',fontsize=18)
plt.xlim(0,80000)
plt.grid(b=True)

plt.tight_layout(1.5,rect=(0,0,1,0.95))
plt.suptitle(u'异常值检测与校正',fontsize=22)
plt.show()
